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We consider the effective model of H = —J\ ^2 Ki ■ > Si ■ Sj — Ji 



(Si • Sj) , describing the 



Mott insulating phase with odd number of spin-1 bosons in optical lattices (J2 > Ji > 0). In terms 
of an SU(3) boson representation, a valence bond mean field theory is developed. In ID, a first- 
order quantum phase transition from a spin singlet to a spin nematic phase with gapful excitations 
is identified at J1/J2 = 0.19833 , while on a 2D square lattice a spin nematic ordered phase with 
gapless excitations prevails. In both ID and 2D cases, we predict that the spin structure factor 
displays dominant antiferromagnetic fluctuations, while the quadrupole structure factor exhibits 
strong ferroquadrupolar correlations. 

PACS numbers: 03.75.Mn, 32.80.Pj, 75.10.Jm 



Degenerate alkali atoms are considered weakly inter- 
acting boson gas due to the smallness of the scattering 
length compared with the inter-particle separation. How- 
ever, the situation changes dramatically when an optical 
potential created by standing laser beams confines parti- 
cles in valleys of the periodic potential and strongly en- 
hances local interactions. Recently, the nontrivial Mott 
insulating state of bosonic atoms in optical lattices has 
been demonstrated experimentally 0. In addition to 
solid state systems, spinor atoms in optical lattices pro- 
vide a novel realization of quantum magnetic systems 
with the possibility to tune various parameters of the 
effective models in the absence of disorder 0, 0, H| . 

Alkali atoms have a nuclear spin 3/2. Lower energy 
hyperfine manifold has three magnetic sublevels and a 
total moment 5=1. In order to observe the quantum 
spin phenomena experimentally, one has to consider a 
system with small number of particles and strong inter- 
actions. In the insulating state, atoms are localized, and 
fluctuations in the particle number on each site are sup- 
pressed. Virtual tunnelling of atoms between neighboring 
sites induces effective spin interactions, leading to novel 
quantum magnetic phases [a, Ifl III l3 ■ 

The boson-Hubbard model is used to describe the low- 
energy physics of spin-1 bosonic atoms in an optical lat- 
tice Km, 



H = -t ( a {m a j,™ + h.c.) + -£-^2ni(ni-l) 



where n % = ~£ m % m a i,m and Sf = J2 m ,n a i,m T m,n a i,n 
denote the number of atoms and the spin operator on the 
i site (m = —1,0,1), respectively, while T a denote the 
matrices for the spin-1 particles. 

In the t = limit with odd number of bosons per 
site, the bosonic symmetry of the wave function requires 



that each site has a localized spin with odd quantum 
numbers. The interactions of Eq.(l) are minimized when 
the localized spins take the smallest possible value S = 
1. Away from this limit, the second order perturbation 
theory in t leads to an effective spin superexchange model 

a 0,1a 



H = — Ji 
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where the coupling parameters are confined to < 
J1/J2 < 1. For a two-spin problem, the first term fa- 
vors a polarized configuration, while the second prefers 
a spin singlet form. Thus, the effective model acts as 
a spin-1 Heisenberg antiferromagnet with ferromagnetic 
frustrations. 

It has been an outstanding issue to study the ground 
states of this effective model H H B |H 0- More ex- 
otic ground states may be realized in ID |8Ul3tll4j. but 
there is no consensus yet. In this paper, using an SU(3) 
boson representation, we develop a valence bond mean 
field (MF) theory. In ID, a first-order quantum phase 
transition is found from a spin singlet (J1/J2 < 0.19833) 
to a (short-range ordered) spin nematic phase (J1/J2 > 
0.19833) with gapful excitations, while on a 2D square 
lattice an ordered spin nematic phase is always obtained. 
Moreover, the spin and quadrupole correlation spectra 
in both disordered and ordered nematic phases are cal- 
culated explicitly, to be checked by future experiments. 

To describe the spin-1 operators, an SU(3) bosonic rep- 



(1) resentation is defined by generators F^(i) 



where indices m and n specify the spin projection with 
values — 1,0, 1. The commutation relation is satisfied 

[KM, FW)] = ki [*wK(i) - . (3) 

which forms an SU(3) Lie algebra. The corresponding 
spin operators are expressed as = v / 2(<zj a i,-i + 
a} ;1 a ii0 ), 57 = (S+y, Sf = (aTjOj,! - aT_ 1 a i ,_i). With 
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this boson representation, we can verify [5+, 5.] — 
2S , ?S i;i , [Sf,Sf] = ±S?6ij. In order to fix Sf = 

S(S + 1) = 2, a local constraint J2 m a im a hm = 1 has 
to be imposed. 

To make the spin representation symmetric, linear 
combinations are introduced: bi.i = A= (a,i,-i — cti,i), 

6»,2 = ( a i,-i + a -t,i)) and & i,3 = Gi.Oj and the spin 
(dipolar) operators are written in antisymmetric form 



Sf 
Sf 
Sf 



-i(bl 2 bi,3 

-*( 6 i,3 6 i,l 
A,2 



&!,3^,2), 
&I,l&i,3), 



(4) 



and the local constraint holds as V foj 6 - = 1. Actu- 
ally for spin-1 bosonic atoms, quadrupole operators can 
be defined by 



? (0) 
,(2) 



= 3(5? f - 2 = 6+6,,! + 6t 2 6 ii2 - 2b{ 3 6 i>3 , 
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(&I,2^,3 + &J, 3 &i,2). (5) 



Three dipole and five quadrupole operators form gener- 
ators of the SU(3) Lie group, as in the Gell-Mann ma- 
trix representation. Then, the effective spin model is 
expressed as: 



<i,j> at, (3 

+ (J 2 - Ji)b{ a hAh 3 ,p ■ 



(6) 



In the limit J\ = J2, the model is reduced to an SU(3) 
symmetric ferromagnetic superexchange model, invari- 
ant under the uniform SU(3) transformation, for which 
the spin ferromagnetic and spin nematic long-range or- 
der can coexist in higher dimensions |l5j . On the other 
hand, in the limit J\ — 0, the model is reduced to an 
SU(3) symmetric valence-bond antiferromagnetic model, 
invariant under the staggered conjugate transformations 
of the two sublattices, leading to a spin dimerization in 
ID p[I3. 

It is known that an SU(2) Schwinger boson MF theory 
can describe rather well the spin correlations for spin-1/2 
Heisenberg antiferromagnets |18|. To develop a similar 
valence bond MF theory, we introduce symmetric pairing 
parameters 
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(7) 



To preserve the SU(2) spin rotational symmetry of the 
model, we have to assume A^i = A 2j 2 = A3 j3 = A a 
and Aj.,2 = A1.3 = A 2j 3 = A&. By introducing a Nambu 
spinor * k = (&£ v 6 k2 , 6 k3 , 6_ k ,i, 6_ kj2 , 6_m), we can 
rewrite the MF Hamiltonian in a matrix form 



+ 1 [(3J a - Ji) A 2 a + JiA%] Nz - ^XN, (8) 

where the local constraint has been implemented by a 
Lagrangian multiplier A, N the total number of lattice 
sites, and z the number of the nearest neighbor sites. 
The MF Hamiltonian matrix is given by 

H m/ (k) = A + 2A a (k)cr x <g> 1+ A b (k)a x <g> M, (9) 

1 1 \ 

1 1 , / is the 3x3 unit matrix, o~ a 
110/ 

(a = x,y,z) are Pauli matrices, A a (k) = |A a (3J 2 — 
Ji)7 k , Ab(k) = zA fc Ji7 k , 7 k = i^] 5 cosk-(5, and 5 the 
nearest neighbor vector. The corresponding Matsubara 
Green function (GF) is thus deduced to 



where M = 



G iw n ) = iuj n a z ® I - H m /(k), 



(10) 



where oj n is the bosonic Matsubara frequency. The poles 
of the GF matrix give rise to the quasiparticle spec- 
tra: ei (k) - A/A 2 -4[A Q (k) + A b (k)] 2 and e 2 (k) = 



A 2 - [2 A Q (k) - Ab(k)] , where the lower band ei(k) 
is singly occupied while the higher band e 2 (k) is doubly 
degenerate. 

From the free energy of the system, the saddle point 
equations at T = OK are derived as 
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z[(3J 2 - Ji)A a + 2JiAf, 

2A Q - Aft 
z [(3 J 2 - Ji)A„ - JiA 6 ] ! 



(11) 



ground state energy is E g = 
Ji)A 2 -(- J1A 2 ]. From the GF ma- 
trix, we can also derive the double-time GFs of the 
boson operators 



^zN [(3J 2 



«&k, 2 |& k)1 » = ((6 k , 3 |6 k) i)) = ((6 k , 3 |6 ki2 )) 
- + A) [4A a (k)Aft(k) + A 2 (k)] 



4m k 2 , 



(12) 



and ((6 k ,i|6 kil )) = (<6 k , 2 |6 k>2 » = «6 k , 3 |6 k 3 )>. Using the 
spectral representation, we can calculate the expectation 
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values at T = OK, 



i ^ r a a 

6iV ^ 



(13) 



and (bj^i 



ei(k) e 2 (k) 

1/3. Therefore 



\2 



we find (Sf) = and ((Sf ) z ) = 2/3 (a = x,y,z), 
exhibiting that both time reversal and spin rotational 
symmetries are well preserved in the present MF theory. 
Then the spin nematic state with a quadrupole moment 

(QV) = (Qf z ) = (QY) = ~Q * s a possible ground state. 

To evaluate the spin spatial correlations in the Mott in- 
sulating phase, the dynamic correlation functions should 
be calculated. By expressing the spin operators in terms 
of the Nambu spinor, the spin correlations are given by 

X a,/3 (q, iw m ) 

= X! Tr I r " G ( k ' * w n)r^G(k + q, iuj n + iw m )] , 

where T x , T y , T z denote the corresponding 6x6 matrices 
of Sf, Sf , Sf, respectively, and both w„ and tu m are 
bosonic Matsubara frequencies. Inserting the GF matrix 
and tracing over matrices, we find 



X X ' X (<1, iUn 



(14) 



displaying an SU(2) spin rotational symmetry. The 
corresponding imaginary part x"(q, ui) can also be ob- 
tained through performing the summation over the 
Matsubara frequency and analytic continuation |l9j| . 
According to the fluctuation-dissipation theorem, the 
static spin structure factor is thus obtained So (q) = 
— duj [1 + jig(w)] x"(q, cj). Furthermore, when the 
correlation function of the quadrupole operator is defined 
by S Q (r-r') = ([3(S r z ) 2 - 2][3(S r z ,) 2 - 2]), the static 
quadrupole structure factor Sq (q) can be evaluated as 

well [la. 

In ID, z = 2 and 7k = cosfc. For a given value 
of Ji/J 2 , the saddle point equations are solved numer- 
ically, results being displayed in Fig.l. For < J1/J2 < 
0.19833, we find A 6 = 0, A a = 0.5077 J 2 , the quasipar- 
ticle band ei(fc) and £2(fc) are degenerate, so the ground 
state is a spin singlet with an energy gap. From the point 
of view of spin correlations, such a spin singlet phase is 
similar to the spin dimerized state in the limit J\ = 
H3. For J1/J2 > 0.19833, both A 6 and A a are finite, 
and two energy gaps in the quasiparticle bands are found 

A s , x = ^/A 2 - 4 [(3J 2 ~.h)A a + 2,hA b } 2 and A Sj2 = 

■\J \ 2 — 4 [(3J2 — Ji)A a — JiA;,] 2 at momenta k = 0, ±7r. 
Moreover, the quadrupole moment Q jumps from to 
0.16897 at the critical coupling Ji/J 2 = 0.19833, and 
then the system is in a short-range ordered spin singlet 
phase with a finite quadrupole moment. Thus, there is a 
first-order quantum phase transition from the spin sin- 
glet to spin nematic phase with gapful excitations. 
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FIG. 1: Saddle point solution in ID is displayed: the valence 
bond parameters (a), the ground state energy (b), two energy 
gaps of quasiparticle bands (c), and the quadrupole moment 
(d) as functions of J1/J2. J2 is taken as the energy unit. 



If we rewrite J\ = — cos 6 and 3 2 = — sin 9, the crit- 
ical coupling corresponds to 6 C = — 0.56337T. In fact, 
a disordered spin nematic (non-dimerized) phase with 
breaking the SU(2) spin rotational symmetr y w as first 
suggested around 9 ~ —3tt/4 by Chubukov [T^]. How- 
ever, subsequent numerical work did not support 
this proposal and it was then believed that the dimerized 
phase prevails in the whole regime — 37r/4 < 9 < —tt/2, 
i.e., < J1/J2 < 1. On the other hand, recent quan- 
tum Monte Carlo l2Cl . density matrix renormalization 
group calculations j2J, and quantum field theory ap- 
proach [22l | indicate that the dimerized phase may end 
at 9 ~ — 0.677T, casting doubts on the earlier anticipa- 
tions. 

To put our MF results on a solid ground, the static 
spin structure factor So (?) and quadrupole structure fac- 
tor Sq(q) are calculated at T = OK and displayed in 
Fig. 2. We find that Soiq) exhibits a broad peak around 
the antiferromagnetic wave vector q = it in both spin 
singlet and gapped spin nematic phase. However, the 
quadrupole structure factor Sq{q) shows a broad peak 
around the ferromagnetic wave vector q = in the spin 
singlet phase, while SQ^q) exhibits a sharp resonance at 
q = in the spin nematic phase, indicating strong ferro- 
quadrupolar spatial correlations. 

On a 2D square lattice, we have z — 4 and 7k = 
(cosk x + cosk y )/2. At T = OK, the conversion from 
summation over momentum to integral will be invalid 
as A — > 4 [(3J2 — Ji)A Q + 2JiAb], and then the quasi- 
particle band ei(k) becomes gapless and linear near 
k* = (0,0) and (7r,7r). The Bose-Einstein condensation 
thus occurs, while a finite energy gap still exists in the 
band e 2 (k). By separating the divergent term from the 
summation, we can introduce a superfluid density p s and 
then solve the saddle point equations. The results are 
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FIG. 2: Spin structure factor (a) and quadrupole structure 
factors (b) in ID for spin singlet phase (solid line J1/J2 = 0.1) 
and spin nematic phase (dashed line J1/J2 = 0.5). 
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FIG. 3: Saddle point solution on a 2D square lattice is shown: 
the valence bond parameters (a), the superfluid density and 
quadrupole moment (b) as functions of J1/J2. 

shown in Fig. 3, where a finite quadruple moment Q is 
always obtained and it increases as the coupling J1/J2 
grows. The ground state in the the regime < J1 / J2 < 1 
is thus a long-range ordered spin nematic phase. 

It is also more important to examine the spin and 
quadrupole spatial correlations in the ordered nematic 
phase. The calculated static spin structure factor 
S jj {q x , q y ) is shown in Fig. 4a, where a strong broad peak 
is displayed at q = (jr, tt), indicating dominant antifer- 
romagnetic spin correlations, and weaker broad peaks 
at q — (0, 7r) and q= (tt,0), corresponding to collinear 



spin-dimer correlations. As a comparison, the calculated 
static quadrupole structure factor Sq (q x , q y ) is displayed 
in Fig. 4b, where a (5-like resonance appears at q = (0, 0), 
implying strong ferroquadrupolar long-range correlations. 
These are unique features of an ordered spin nematic 
phase, which can be directly probed by polarized inelas- 
tic light scattering experimentally. 

To summarize, the SU(3) boson representation is used 
to develop an efficient valence bond MF theory for the 
Mott insulation phase with odd number of spin-1 bosons 
in optical lattices. A first-order quantum phase transi- 
tion from a spin singlet to a nematic phase with gap- 
ful excitations is predicted for ID, while on a 2D square 
lattice a spin nematic ordered phase is shown to always 
prevail. Both predictions are further supported by an ex- 
plicit calculation of the spin and quadrupole correlation 
functions. 
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